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The one-electron spectral-weight function of the half-filled ionic Hubbard model 
is calculated by means of Quantum Monte Carlo. A metallic regime occurs between 
two values of the coupling constant (Hubbard U) U\ < XJ%. The system is a band 
insulator below Ui, with a strong charge density wave corresponding to a large 
ionicity, and an increasingly antiferromagnetic (AF) insulator above U2 evolving 
into a Mott insulator as U — > 00 . The intermediate regime, which is both AF and 
dimerized, is caused by a two-peak structure at A;a=0.5-7r (the Fermi surface). As U 
increases both peaks approach each other, overlap, and separate again, the system 
becoming metallic in the overlap region. This behavior can be traced to a site-parity 
change of the ground state(GS), corresponding to a curve-crossing between two GS 
of different site parity. 

PACS numbers: 71.30,+h,71.10. Pm, 77.80.-e 

I. INTRODUCTION 

In quasi-one dimensional materials like halogen-bridged transition metal chain complexes, 
conjugated polymers, organic charge transfer salts, or inorganic blue bronzes, the itinerancy 
of the electrons strongly competes with electron-electron and electron-phonon interactions 
which tend to localize the charge carriers by establishing commensurate spin or charge- 
density wave ground states (GS). At half filling band (BI) or Mott (MI) insulating phases 
are favored over the metallic state. Quantum phase transitions between the insulating 
phases are possible and the character of the electron excitation spectra reflect the properties 
of the different GS. A controversial issue is the nature of the BI-MI transition as well as 
whether or not only one critical point separates both insulating phases in purely electronic 
model HamiltoniansBi. Despite their importance in low dimensional materials!'!^ phonon 
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dynamical effects will be ignored in this paper which will be concerned with the so-called 
one dimensional (ID) ionic Hubbard model (IHM), just a Hubbard model with a staggered 
potential which displaces by A the energy levels of even and odd sites. It is described by 
the Hamiltonian 

<ij>s is i ^ ' ^ * 

where electrons of spin s — ±| hop with amplitude t between nearest-neighbor inequivalent 
sites of energy level ±|A(A = lin what follows). Electrons of opposite spin on the same 
site feel a repulsion U. 

This model Hamiltonian, originally proposed by Nagaosai and later by Egamii as a 
model ferroelectric, is ideal for studying the issue of quantum phase transitions in electron 
systems. On general grounds one expects a transition to occur from an ionic band insulator 
to a strongly-correlated Mott insulator as U increases. Evidence for such a transition was 
found in exact-diagonalization calculations by Resta et ali predicting a metallic point at a 
critical value U = U c separating both insulating phases. This transition is signaled by a 
polarization jump with a sign change associated with a site parity change of the GS. The 
same conclusion was supported by Gidopoulos et afl, who showed the reversal of site parity 
above to be of magnetic origin, thus leading to a vanishing of the spin gap, A s = 0. Likewise, 
Brune et ali by means of the bosonization technique along with Lanczos and density matrix 
renormalization group (DMRG) calculations lent support to this picture. 

On the other hand, Fabrizio et ali shed new interest into the IHM in a field-theoretical 
bosonization analysis where they propose a new scenario for the intermediate region, U ~ A, 
between the BI and the MI. The model should instead exhibit two quantum phase transitions: 
one from a BI state to a long-range bond-ordered (BO) state predicted to be in the Ising 
universality class and a second one from the BO to the MI state predicted to be a Kosterlitz- 
Thouless transition. Such transitions to BO states have recently been found in ID Hubbard 
models with extended interactions (U-V) by Nakamurai. The BO state is a broken symmetry 
state in which the system becomes ferroelectric due strictly to electron-electron interactions 
even if all the atoms are at centers of inversion. This is also called the spontaneously 
dimerized insulating (SDI) phase, which, as argued by Fabrizio et ali, implies a finite spin 
gap, As >- 0. Such a SDI phase has been observed numerically by Wilkens and MartinS 
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in a variational quantum Monte Carlo study. More precisely, they found one transition 
from the BI to the correlated SDI phase, but no second transition to the MI phase which 
is approached asymptotically by the SDI phase as U — > oo. Recent DMRG calculations0ll! 
lead to mutually conflicting results. Thus whereas Takada and KidoS support the Resta et 
al@ picture of a single first-order transition, Qin et al0 support the Fabrizio et afl scenario 
of two second-order quantum phase transitions. 

In this paper we present the one-electron spectral-weight function A(kcj) of the IHM. 
The calculation is done by means of a Quantum Monte Carlo (QMC) simulation in the 
grand canonical ensemble suplemented by an approximate effective action. Its eigenvalues 
and eigenvectors alow us to construct explicitly dynamic correlation functions in real time or 
frequency. This procedure leads to spectral functions with a reasonably high resolution. The 
calculation is performed at a relatively large inverse temperature, [3t — 10 (t — 1 in what 
follows), which brings the system very close to its GS properties. Our calculated A(kci;)clearly 
shows metallic character (non- vanishing density of states in the gap) at ka = 0.5tt (a = 1), 
between two critical U values instead of a single metallic point at the first threshold. Thus a 
double BI-M-MI smooth transition is predicted, the intermediate metallic regime at U~ 2A 
being both antiferromagnetic (AF) and dimerized 

II. EFFECTIVE HAMILTONIAN AND CORRELATION FUNCTIONS 

In our calculations, we use the QMC approach in the grand canonical ensemble, which has 
been explained in detail by Hirsch§0 and White et aM As is well-known, the partition 
function is factorized into L time slices of extent At = (3/L. For Hubbard models, the 
interaction U ^ ^t n *i i s transformed into a discrete Ising fielcM an = ±1, which depends 

on lattice site i and imaginary time r = I At (1 = 1, L). The non-interacting Hamiltonian 

for a given Ising configuration (a definite allocation of +1 and —1 values for all the an is) 



rea> 



dsS 




(2) 



<ij>s 



where 




(3) 
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The first term on the rhs of Eq(3) is just the hopping term of the IHM, Eq(l). In the second 
term, \i is the chemical potential, a — ±1 for up/down electron spin, and A is a non- linear 
function of U given by cosh A At = exp (C/Ar/2).In the limit At — > the partition function 
can be written as 

Z=Y J TrYi l e- A ^ = J2 Tr e ~ ATEl H ' = J2 Tr z- m+V) (4) 

<cr> <cr> <<T> 

where J2<a> means summing over the Ising configurations and Tr denotes taking the trace 
over the electron degrees of freedom. K is the r-independent part of H(r) and 

V = j a n ( n *T - n u) = A ^ < ^ > (n n - nn) (5) 

il i 

< <Ji >= (1/L) a u being the average Ising spin at the i-th site. Therefore, H e ff = K + V 
is r-independent in the limit At — > 0. We adopt this effective r-independent Hamiltonian, 
for doing measurements even if the number of time slices used in the calculation is finite. 
To minimize as much as possible size effects due to a finite L, we shall increase L and seek 
for convergence. Needless to say, the update of the Ising field is done following the standard 
procedure of, say, Ref 15. 

The corresponding Hamiltonian matrix for each spin s (Eq(3) with < a >instead of 

< an >) is now diagonalized and, in terms of its eigenvalues e MS and eigenvectors | /is >, we 
can construct correlation functions in real time or frecuency. For instance 

for the causal Green's functions and likewise for other quantities. As customary, f^ s = 
(e 136 ^ + l) 1 is the Fermi-Dirac distribution and rj a vanishingly small imaginary part. All 
these quantities must now be averaged over the accepted Ising configurations to obtain the 
final result. It is just in this averaging process where the correlation effects are restored. 
Since the final result is translationally invariant, one can Fourier transform and find the 
spectral-weight function (SWF) 

A ks {^) = --Im{G>M ~ G< s {u)} , (7) 

71 

and the DOS 

k> 
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which can also be obtained from (6) as (^/N)'^2 i N i (uj), N^co) being the corresponding 
imaginary part of G s u (uj). For our case of the IHM, this last route allows to calculate the 
partial DOS for even/odd atoms by simply restricting the sum over even/odd sites. As 
noticed earlier, this building procedure leads to a good resolution of spectral features. 

III. DOUBLE INSULATOR-METAL-INSULATOR TRANSITION 

We characterize the system by its one-electron SWF, A k (uj), along with the charge and 
spin structure factors, c k and s k , defined as the Fourier transform of the static charge and 
spin, respectively, correlation functions. These are given as usual by =< q^qj > and 
Sij =< s iz Sj Z >, qi = n^+riii and s iz = (l/2)(n^— n^) being the charge and spin at the i—th 
site and r = 0. This type of characterization is useful and can be of help in understanding 
the electronic properties of a system. For instance, the SWF tells us the presence or absence 
of gaps or pseudogaps as U increases, which is of fundamental importance in establishing 
the insulating or metallic character and possible transitions between them of, say, the IHM. 
Likewise, the evolution of c k and s k with U, especially theirs peak at k — ir, gives information 
about the charge or magnetic origin of the transitions. 

A. The spectral-weight function. 

As hinted before, it is just the structure of A k (cu) at k = it/2 which accounts for the 
insulating or metallic character of the IHM. Hence, Fig 1 shows this quantity for three 
values of U (recall that t — 1) so selected to exhibit the three fundamental regimes of the 
model. The calculation has been done for a chain of eighty sites with periodic boundary 
conditions at half filling and f3 = 10. Typically 150 time slices, two hundred warm sweeps 
and one thousand measurements were taken. At U — 1, a two-peak structure with a clear 
gap in between signals an insulating system. At U — 2.25, both peaks have overlapped 
developing a third peak just at the midpoint of the gap (which becomes the Fermi level). 
This three-peak structure clearly shows a metallic system. Finally, at U — 4 the central 
peak splits into two peaks leaving a gap in between. This four-peak structure signals again 
an insulating system which must be somehow of a different kind from the initial one (see 
below). As U increases (not shown) the central gap widens and the four peaks go over into 
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FIG. 1: The spectral-weight function A^{u) for k = ir/2 and U/t (t=l) values lying in the three 
coupling regimes of the ionic Hubbard model. The QMC calculation has been made for a periodic 
chain of eighty sites an j3 = 10 



a two-peak structure as U — ► oo, the system becoming then a Mott insulator. 

In order to estimate the approximate extension of the metallic region, Fig 2 shows the 
evolution of the spectral density at k = n/2 and oo = (the midgap) as U increases. We see 
a maximum at U m = 2.25 an onset threshold at U\ ~ 1.75 and a falloff with an approximate 
endpoint at U 2 = 2.75 eV leaving a zone of metallicity U2 — U1 = 1. Since this quantity may 
be size-dependent, we proceed to a finite-size scaling in order to approach the thermodynamic 
limit. Calculations made on chains of 8,12,16,80,120, and 160 sites show that U2 — U\ — 1 is 
roughly independent of size, but U m tends to go up, extrapolating to 2.50 when 1/N — > 0. 
Notice that only finite chains of iV = Am sites should be used since only them have k = ir/2 
as an allowed wave vector. Chains of other lengths fail to detect the metallicity. In the 
thermodynamic limit, of course, all chains have a wawe vector innnitessimally close to ir/2 

The three-peak structure responsible of the metallic phase disappears as soon as one 
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FIG. 2: Ak{oj) at k = it/2 (a=l) and E=0 (midgap which becomes the Fermi level) as U increases. 
Same calculations of Fig. 1 



departs from k = ir/2, showing a single mainly occupied peak at k = (7r/2) — 5 and a mainly 
empty peak at k — (ir/2) + 5 ( 5 = 2it/N, N being the number of sites). Hence the metallic 
peak should be visible in angle-resolved photoelectron spectroscopy and its inverse, but not 
in integrated measurements which give the DOS, Eq (8). Since the weight of that structure 
falls with iV (Fig 3), only a pseudogap with a sizeable DOS at the Fermi levels should be 
seen. 



B. The charge and spin structure factors 

Fig 4 finally displays the evolution of the peak of and s fc at k = it for increasing U. 
The charge factor consistently drops from a high value at U — 1, which signals a strong 
CDW in correspondence with a large ionicity (the fractional ionic character of the bond is 
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FIG. 3: Density of states for the same chain as Fig 1 for U and number of sites N=8,12,16 and 80. 

FIC = 0.32) which also falls down monotonically. The system, therefore starts as an ionic 
band insulator at small U. The spin structure factor starts close to zero, corresponding to 
a paramagnetic insulator, and monotonically increases with U, signaling an AF SDW of 
increasing amplitude. For large U, the systems evolves into a Mott insulator with a strong 
AF SDW and a vanishingly small CDW, i.e., zero ionicity (neutral). 

In the intermediate regime, U\ < U2, the system is both increasingly AF and neutral 
(decreasing ionicity, i.e., decreasing CDW). On the other hand, Table I displays a small 
sector of the density matrix (nfj =< cf s Cj S >), for U = 2.25 just the top and left 4x4 
matrix, which is sufficient to show that the metallic state of this intermediate regime is 
dimerized. Thus n s 12 > n 23> e ^ c > indicating that the sites (12), (34), and so on, form dimers. 
These dimers weaken for U > U 2 and tends to disappear as U — > 00 
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FIG. 4: The peak at k = it of the charge, (triangles), and the spin, Sk (circles), structure factor, 
as U increases, for the 80-sites periodic chain. 

IV. CONCLUDING REMARKS 

Our QMC study shows three coupling regimes in the ionic Hubbard model: (i) An ionic 
band insulator below U = U\.(ii) An increasingly antiferromagnetic and dimerized metal 
which is gradually becoming neutral in the intermediate coupling regime Lq < U < U 2 - And 
(in) an increasingly AF and neutral insulator whose dimerization, measured by n 12 — n 23 , is 
gradually dropping away for U > U 2 . This last regime tends to a Mott insulator a U — > oo. 

A remark is here in order. Although the spectral- weight function at k = n/2 seems to 
hint a critical separation among the three coupling regimes (a threshold and an endpoint), 
the smooth behavior of all the integrated quantities (DOS, ionicity, energy, charge and spin 
structure factors, etc.) rather suggest a smooth interpolation between the three regimes 
which should be distinguished from true quantum phases with sharp, critical, separation 
among them. 
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To conclude, let us finally say that the results of our QMC calculations would be quite 
similar to those of Resta et alS, Wilkens and Martin0,etc, if the metallic zone reduced to a 
single point. On the other hand, they would also be quite close to those of Fabrizio et all 
should the intermediate regime reported in this letter turn out to be insulating. 

We acknowledge the financial support of the Spanish DGICYT through Project NPB 
98-0683 
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TABLE I. Upper-left 4x4 sector of the density matrix nf ■ =< c^c,- s > 



0.5477 0.2928 0.0109 -0.0555 

0.2928 0.4523 0.2116 -0.0133 

0.0109 0.2116 0.5477 0.2928 

0.0555 -0.0133 0.2928 0.4523 



